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ABSTRACT 

We present simulations of the tidal disruption of a solar mass star by a 10 6 Mq black hole. These, for the first time, 
cover the full time evolution of the tidal disruption event, starting well before the initial encounter and continuing 
until more than 90% of the bound material has returned to the vicinity of the hole. Our results are compared to the 
analytical prediction for the rate at which tidally-stripped gas falls back. We find that, for our chosen parameters, 
the overall scaling of the fallback rate, Mfb, closely follows the canonical t -5 / 3 power-law. However, our simulations 
also show that the self-gravity of the tidal stream, which dominates the tidal gravity of the hole at large distances, 
causes some of the debris to recollapse into bound fragments before returning to the hole. This causes Mfb to vary 
significantly around the t -5 / 3 average. We discuss the implications of our findings in the context of the event Swift 
J1644+57. 

Subject headings: black hole physics — galaxies: nuclei — X-rays: individual (Swift J1644+57) — 
hydrodynamics 


1. INTRODUCTION 

When a star comes within a supermassive black hole’s 
(SMBH) tidal radius r t ~ iL^M^/M*) 1 / 3 , where R* 
and M* are the stellar radius and mass, respectively, 
and Mh is the black hole mass, the tidal force exerted 
by the hole across the star is sufficient to overcome its 
self-gravity, resulting in its destruction. Early studies of 
these tidal disruption events (TDEs) demonstrated that 
roughly half of the disrupted ste l lar debris is bound to 
the black hole (jLacv et al.l 119821: lR,eesl 1 19881 ) , meaning 
that it will eventually return to the tidal radius, dis¬ 
sipate energy through shocks (lEvans fe Kochanekl 119891: 
iKochanekl 11994 IGuillochon et al.l l2014bl ). and form an 
accretion disk. The resulting accretion power is then ca¬ 
pable of producing a highly lu minous event, and many 
have already been detected (iKo mossa & Bade 1999; 
Gezari et al.l 120081: iBurrows et al.l 1201 It iCenko et al.l 
20121: iBogdanovic et al1l2014l). 

Considerations of the star at the time of disruption 
show that the rate at which the bound, tidally -stripped 
mate rial returns to the hole decreases as £ -5 / 3 (jPhinnevI 
H989lh and early simulat ions supported this estimate 
(jEvans fc Kochanekl 11 9891 ). More recently, authors in¬ 
vestigated the consequ ences of the stellar composition 
on the rate of return (jLodato et al.ll2009 h illustrating 
that the early fallback stages have a more complex tem¬ 
poral behavior that depend on the density stratification. 
IGuillochon fc Ramirez-Ruiz ( 2013 ) showed that the peri- 
center distance can have dramatic effects on the rate 
of return of debris. General relat ivistic Lense-Thirring 
(jStone &; Loebll20l2l ) and apsidal (IBonnerot et al.ll2015l) 
preces sion and recompression shocks (jGuillochon et al.l 
l2014bl ) have also been investigated. 

Most of the recent work on TDEs has been per¬ 
formed with__thu_^^ simulations. For ex¬ 

ample, IGuillochon et al.l (j2014b[ ) used a grid code (see 
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also Shiokawa et al. 2015) to follow the evolution of the 
tidally disrupted debris from the initial encounter with 
the black hole to its e ventual return to pericenter, while 
iBonnerot et al.l (j2015f ) used a smoothed particle hydro- 
dy namics (SPH ) code to achieve the same feat (see also 
lHavasaki et al.l 20151b In these cases, however, the set 
of parameters used to model the tidal disrupt ion was 
somewhat unphysical - IGuillochon et ahl (|2014bl ) used a 
10 3 Mq black hole when followin g the return of the de- 
bris to the pericenter radius, while IBonnerot et al.l (j2015f ) 
used a more tightly bound stellar progenitor (e = 0.95 
for their most parabolic simulation, e being the stellar 
eccentricity). Both of these choices were made to reduce 
the computational cost at the expense of reality. Instead, 
a more likely mass for an SMBH is 10 5 - 10 8 Mq, and 
the infalling star is initially so far from the hole that its 
orbit is effectively parabolic (e ~ 1). 

A full TDE with realistic parameters has not yet 
been modeled because of the extreme set of spatial and 
temporal scales. In spite of these difficulties, we present 
the first simulation that has resolved the full duration of 
a TDE with a solar progenitor on a parabolic orbit, the 
periapsis of which coincides with the tidal radius, and 
a 10 6 Mq SMBH. The simulation runs from well before 
the initial tidal interaction to long after the most bound 
debris returns to pericenter, corresponding to roughly 
ten years after disruption. This enables us to calculate 
explicitly the rate of return of tidally stripped debris to 
the black hole, demonstrating the effects of self-gravity 
which dominate that of the hole when the material is 
near apoapsis. In section 2 we provide analytic formulae 
for the shape and density distribution of the stream as it 
recedes from the hole, and we develop arguments which 
suggest the dominance of self-gravity over the tidal 
field of the hole at large radii. In section 3 we describe 
the initial conditions and numerical method used to 
simulate the encounter. Section 4 presents the results 
of the simulations and their analysis. We conclude and 
discuss further implications of our findings in section 5. 
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2. SELF-GRAVITY 


The forces acting on the debris are self-gravity, the 
tidal shear across the width of the stream, and pressure. 
If we consider a small volume within the stream of length 
SR and mass SM , then the self-gravity of the stream 
dominates the tidal field if 


Po 


3_Mh 
27r r 3 


(i) 


where po = 3SM/ (MSR 3 ). 

We can calculate when this condition is satisfied by 
approximating the stream as a circular cylinder of width 
H. In this case, the differential amount of mass dM 
contained in the stream is 


dM = irH 2 pdx , ( 2 ) 

where x is the length measured along its center. If the 
specific energy distribution of the gas is frozen in at dis¬ 
ruption, then the curve traced out by x is a function 
solely of the position of a given gas parcel when the star 
is at pericenter, which we will denote by R p = /uR*, 
with — 1 < p < 1. Modeling the star as a polytrope 
([Hansen et al.ll2004L the di f ferential amount of mass con¬ 
tained in dp, is ([Coughlin fe Begelmanl 120141 ) 


m i j^em# 

^ 2 J^ 1 0 (O n £ 2 d£’ 


(3) 


where n is the polytropic index of the gas, 0(0 is the 
solution to the Lane-Emden equation and 0 is the first 
root of 0 ( 0 - 

We can parametrize x by the functions r(p) and 0(p ), 
where r and 0 are solutions to the equations of motion 
of a test particle in the Newtonian gravity of the hole: 

£ 2 1 () 
GMh 1 + e cos 6 ’ 

r 2 0 = £. (5) 

Here £ and e are the specific angular momentum and 
eccentricity, respectively, of the gas parcel with position 
//, and are given by 


i = v / 2GM/ l r t (l - nq 1/3 ), ( 6 ) 

e = l-2 (7) 

where q = We have assumed here that the gas 

parcels initially move with the center of mass of the star 
and /3 = r t /r p = 1. By numerical integration of equation 
© we can determine 0(p ), and r(p) is then found by 
using that result in equation 

T he infinitesimal distance along the curve is dx = 
dp^J (r ') 2 + r 2 (#') 2 , where primes now denote differen¬ 
tiation with respect to p. Using this result and equation 
© in equation © gives 

M i & (a 

27rffVM 2 +r 2 ( 0') 2 fo 1 e(O n £ 2 d£‘ 

Using this equation in condition © gives the range in p 
over which self-gravity dominates the tidal shear. 


The pressure of the material can also resist self-gravity. 
However, even if a small overdensity of radius SR and 
mass SM is pressure- supporte d, it must satisfy the Jeans 
condition SR < c s yJn/{pG). If SR exceeds this value, 
any perturbation will cause material to clump around 
SM. 

Finally, the local velocity gradient can potentially sta¬ 
bilize the stream. The divergence of the velocity de¬ 
creases the density within the stream on a timescale 
T d iv ^ 1/(V • v). If T div > Tff, where r ff ~ 
SR 3 ' 2 /\/GSM is the local free-fall time, self-gravity will 
overcome the stabilizing effect of the divergence of the 
flow. 

By balancing the tidal shear and s elf-gravity in the 
transverse direction, iKochaneki (Il994h found that H oc 
r 1 / 4 for ay = 5/3 equation of state. We will use this 
when considering how the density (equation [ 8 ]) varies 
with r (see Figure©. 

3. SIMULATION SETUP 

To test our analytical reasoning above and the canoni¬ 
cal £ -5 / 3 fallback rate, we now employ three dimensional 
hydro dynamic simulations. We use the SPH code pha n- 
tom ([Price &; Federrathl 120101 : iLodato fe Pricel l2010h to 
simulate the tidal encounter, phantom was designed for 
running high resoluti on hydrodynamic s im ulations with 
compl ex geometries flNixon et al.l l2012a.b: Ma rtin et al.l 
l2014allbl : iNealon et al.ll2015h . and is therefore well suited 
to studying TDEs. 

For our simulations we model the black hole as a New¬ 
tonian sink particle situated at the origin. Self-gravity of 
the stellar material is modeled using a k-D tree alongside 
an opening angle criterion, ca lculating d irectly the forces 
arising over short distances ([Barnes &; Hutlll986h . The 
opening angle for the simulations presented here was 0.5 
(we ran an additional simulation with a more accurate 
opening angle of 0.3 and found negligible differences with 
the run presented here). The gas retains a polytropic 
equation of state, i.e., only adiabatic energy gains and 
losses are included (we also performed simulations with 
shock heating and, with our parameters , found negligible 
differences; see also lLodato et al.ll2009 ). 

We initialized our runs by placing 10 6 particles in the 
configuration of a 7 = 5/3 polytrope; this is achieved 
by first placing the particles on a close-packed sphere. 
This configuration is then stretched to achieve the cor¬ 
rect polytropic density distribution, the resulting poly¬ 
trope then placed at a distance of 10 r t from the hole 
with its center of mass on a parabolic orbit, its location 
at periapsis equal to the tidal radius (equations (© and 
© with £ = ^/2GMhVt and e = 1). The distance of 
10 rt allows the initial configuration time to relax before 
interacting with the black hole. 

Here we are primarily interested in exploring the ef¬ 
fects of self-gravity on the stream evolution, and so we 
include self-gravity at all stages of the simulation. We 
are also interested in comparing the theoretical rate of 
return of bound material (not necessarily equal to the ac¬ 
cretion rate onto the black hole), Mfb oc t -5 / 3 , with that 
inferred from the simulation. Therefore, once the dis¬ 
rupted debris is beyond periapsis, the accretion radius 
of the black hole, initially well inside the tidal radius, 
is extended to r acc ~ 3 r t (— 120 R g ). Once the bound 
material returns to the point of disruption, it is removed 
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from the simulation. 

4. RESULTS 

The top, left panel of Figure [1] shows the density of 
the debris stream at periapsis, while the top, right is 
0.84 hours af t er dis ruption; these reproduce Figure 4 of 
iLodato et al.l (|2009[ ). The bottom four show, from left to 
right and top to bottom, 5.38 days, 15 days, 22.4 days 
and one month after disruption. The red curves indicate 
the analytic predictions; the fact that they provide good 
fits to the data means that the particles approximately 
trace out Keplerian orbits in the potential of the black 
hole. 

The black points in Figure [2] show the density distribu¬ 
tion along the stream at t = 22.4 days after disruption; 
the solid, red curve gives the analytic prediction (equa¬ 
tion © with H = Rq(v /tt) 1 / 4 ) and the dashed, red 
curve is the density at which the self-gravity of the debris 
and the tidal field of the hole are equal along the stream 
(the right-hand side of equation (pQ); self-gravity domi¬ 
nates for points above this curve). The wings present 
in the si mulated values, wh ich were also found by other 
authors (jLodato et alJl2QQ9h , are due to the fact that the 
shape of the star at the time of disruption is not spherical 
(see the left-most panel of Figure [1]). 

Figure [3] shows the rate of return of material 
to pericenter, Mfb(t), from the simulation (solid, 
black curve) and the analytic solution (red, dashed 
curve). The analytic estimate was determined 
by using equation © and the fac t that fi(t) = 
t Z D- 2 / 5 , where T = 27rM, l /(M*\/GA4)(i?*/2) 3 / 2 

Coughlin fo Begelmanll20T~4t) . The time at which mate- 
rial returns to pericenter is slightly earlier than that pre¬ 
dicted analytically, a consequence of the tidal distortion 
at the time of disruption (see Figure [3]). The asymptotic 
scaling Mjb oc £ -5 / 3 is closely followed by the numerical 
solution at late times, reaffirming that the gas parcels 
follow approximately Keplerian motion. In our simula¬ 
tion, half of the bound material is accreted by about one 
year - shortly after the peak in the fallback rate. 

Figure [3] also shows that the fallback rate quasi- 
periodically deviates from the t -5 / 3 law at late times. 
This scatter is due to the fact that the stream fragments, 
forming small gravitationally bound clumps that enter 
the accretion radius at discrete times, and is not sim¬ 
ply numerical noise. To emphasize this point, Figure 0 ] 
shows the stream two years after disruption and insets 
that focus on clumps near the center of the stream. 

5. DISCUSSION AND CONCLUSIONS 

We have presented the results of a TDE in which a solar 
mass star (with its pericenter at the tidal radius and its 
center of mass on a parabolic orbit) was disrupted by a 
10 6 Mq black hole. Contrary to past investigations, we 
resolved the full duration of the TDE - from the initial 
encounter between the undisturbed star and the hole to 
long after the most bound, tidally-stripped debris has 
returned to periapsis. When finished, more than 90% of 
the bound material was accreted. Crucially, we are the 
first to follow this long-term evolution for a solar mass 
star on a parabolic orbit around an SMBH. 

The two main aims of this paper are (1) to determine 
the effects of self-gravity on the stream far from the black 
hole (see Figures [3] and 0]), and (2) to compare the rate of 


return of bound material to the analytic estimates (Fig¬ 
ure n. The fiducial rate oc t 5 / 3 was determined 
by assuming that the specific energy distribution of the 
material was frozen in at the time of disruption. How¬ 
ever, the self-gravity of the material, which is important 
at large distances from the hole (Section 2), has the po¬ 
tential to alter the specific e nergy distribution (see also 
iGuillochon fc Ramirez-Rui3120131 ). It is therefore nec¬ 
essary to resolve the full disruption process with self¬ 
gravity included at every step to determine the true rate 
of return. 

Our simulations demonstrate that the fallback rate of 
debris closely mimics the theoretically-predicted one, the 
biggest discrepancy arising from the return time of the 
most tightly bound material. Therefore, the material, in 
agreement with intuition, follows approximately Keple¬ 
rian orbits about the hole. Interestingly, however, we also 
found that the fallback rate at late times tends to over 
and under-estimate the t -5 / 3 law, a feature that is due 
to the fragmentation of the stream into gravitationally- 
bound clumps. When one of these clumps is accreted, the 
fallback rate spikes above its average value, while the rate 
is reduced when the lower density debris between clumps 
returns to periapsis. 

We have shown analytically and numerically (see equa¬ 
tion © and Figure [2|) that the local self-gravity of the 
stream dominates the tidal field of the hole at large 
radii. We also find that the pressure within the sphere 
of the center inset of Figure 0 | is p ~ 10 3 dyn/cm 2 . 
Conversely, the pressure necessary for maintaining hy¬ 
drostatic equilibrium is p eq ~ GMp/R ~ AirGp 2 R 2 /3, 
where the radius of the clump in Figure 0 ] is R ~ 0.3 r t 
(= 30 Rq) and the density is p ~ 10 -6 g/cm 3 , which 
gives p eq 10 6 dyn/cm 2 . Finally, the computed diver¬ 
gence of the velocity just outside the clump gives a di¬ 
vergence timescale of r^iv — 1/(V • v) ~ 10 7 s, while the 
infall timescale is Tff ~ 1 /y/Gp ^ 4 x 10 6 s. This shows 
that gas pressure and the local shear of the stream are 
also incapable of overcoming self-gravity. 

The tidal field of the hole, gas pressure, and the local 
velocity shear are all incapable of supporting the clump 
in Figure 0 ] against its own self-gravity. This finding sug¬ 
gests that the original stream was gravitationally unsta¬ 
ble, and small perturbations resulted in its fragmenta¬ 
tion. For the results presented here, these perturbations 
are induced by the discreteness of the numerical method. 
To substantiate this claim, we ran the simulation with 
10 4 and 10 5 particles and found that clumps formed 
sooner with decreasing particle number. However, our 
pressure smoothing length, equal to the gravitational 
smoothing length, was always at least as small as the 
Jeans radius Rj = c s yj tt / (Gp) (se e Figure [5]) ; the col¬ 
lapse was therefore alw ays resolved (jTruelove et al.lH997t 
iBate &; Bnrke rt 1997). Additionally, the first clumps in 
the simulation presented here form around a month af¬ 
ter disruption; at this time, the smoothing length at the 
center of the stream is roughly 0.1 times the width of 
the stream, meaning that the forces are clearly resolved 
at this time (see also Figure 5, which shows that the 
collapse itself is also resolved). Finally, we note that our 
SPH si mulations are not su sceptible to artificial fragmen¬ 
tation (jHubber et al.ll2006 ). and our results demonstrate 
that any physical perturbation imposed on a stream that 
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Figure 1. The star at the time of disruption (top, left), where the colors indicate the column density and distances are measured in units 
of tidal radii (~ 7 X 10 12 cm = lOO it^), 0.835 hours a fter periapsis (top, right; the top two figures coincide with the bottom-left panels of 
Figure 4 a) and b), respectively, of ILodato et al.ll2009l .h and 5.38 days (middle, left), 15 days (middle, right), 22.4 days (bottom, left) and 
1 month (bottom, right) after disruption. The red, dashed curves on the bottom two rows indicate the analytic predictions. 
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Figure 2. The computed density of the stream (solid, black 
curve), the analytic prediction (red, solid curve; equation [S} , and 
the density at which self-gravity becomes important (red, dashed 
curve; equation [lj as functions of r (spherical distance from the 
hole) at 22.4 days from the point of disruption. Points above the 
red, dashed curve are self-gravitating. The numerical solution ex¬ 
tends to smaller radii because the star at the time of disruption is 
not spherical (see also lLodato et alll2009i b 



Figure 3. The fallback rate in solar masses per year calculated 
from the simulation (black, solid curve) and the analytic prediction 
(red, dashed curve) as functions of time in years. The temporal 
shift between the numerical solution and the analytic prediction is 
from the stellar distortion at the time of disruption (see text for 
further details). The variability at late times is due to the accre¬ 
tion of clumps formed in the self-gravitating stream. The variabil¬ 
ity starts about 2 years after disruption; however, this timescale 
depends upon the details of the TDE and its environment, and 
therefore may be shorter. 

satisfies O is unstable to collapse. 

To further investigate the instability of the stream, we 
performed additional runs with 10 5 and 10 6 particles and 
seeded the initial distribution of debris with a small, but 
resolved, perturbation. In these cases, clumps formed at 
almost the exact same time with nearly identical masses. 
Therefore, a resolved noise field imposed on top of the 
otherwise smooth distribution of stellar debris results in 
converged fragmentation. This finding further supports 
the fact that the stream is gravitationally unstable, and 
any small physical perturbation to the distribution of the 
debris will cause it to fragment. We will present a more 
complete analytical and numerical analysis of the gravi¬ 



Figure 4. The density two years after disruption and, in the 
insets, close-ups of the fragments that forms along the stream. 
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Figure 5. The ratio h/Rj , h being the smoothing length and Rj 
the Jeans radius, as a function of density at a time when the first 
clumps start to form. 


tational stability of the stream in a future publication. 

Real, fragmentation-inducing perturbations could be 
caused by a number of physical processes. For example, 
the interaction of the debris with any material surround¬ 
ing the hole, the density distribution of which is neither 
necessarily smooth nor homogeneous, could cause local 
deformities within the stream. In light of this, the vari¬ 
ability in the light curve of the event Swift J1644+57 
(iBurrows et al.l l20TTt iLevan et~ahl 120111: IZauderer et al.l 
120111) . the recently-observed, jetted TDE, could be in¬ 
terpreted as abrupt changes to the accretion rate in¬ 
duced by the fallback of bound clumps. This notion 
is supported by the lack of any ultraviolet and optical 
emission from the event, which is indicative of a large 
amount of dust presen t in the circumnuclear environment 
([Burrows et al.ll201lD : interactions between the stream 
and this natal environment could have initiated early 
fragmentation of the stream, causing variability in the 
fallback rate on a timescale commensurate with observa¬ 
tions. 

Finally, we note that while the mass contained in the 
fragments in Figure|4]is SM ~ 0.005Mq, there is a range 
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of clump masses. Therefore, we suggest that the ob¬ 
ject G2 - t he clump of materia l observed near the galac¬ 
tic center (jBurkert et al.112012 ) - and other such clouds, 
could have been pro duced by the tidal disrupti on of a star 
in the recent past. iGuillochon et al.l (j2Q14al ) reached a 
similar conclusion; their clump, however, was formed by 
fluid instabilities generated through the interaction of a 
stream of debris with an ambient medium. Here, on the 
other hand, self-gravity causes the stream to fragment. 
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